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ABSTRACT 

A method for accurately solving inviscid compressible flow in the subcritical and supercritical 
regimes about complex configurations is presented. The method is based on the use of unstruc- 
tured triangular meshes in two dimensions, and special emphasis is placed on the accuracy and 
efficiency of the solutions. High accuracy is achieved by careful scaling of the artificial dissipa- 
tion terms, and by reformulating the inner and outer boundary conditions for both the convec- 
tive and dissipative operators. An adaptive grid refinement strategy is presented which 
enhances the solution accuracy for complex flows. When coupled with an unstructured mul- 
tigrid algorithm, this method is shown to produce an efficient solver for flows about arbitrary 
configurations. 
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No. NAS1-18107 while the author was in residence at the Institute for Computer Applications in Science 
and Engineering (ICASE), NASA Langley Research Center, Hampton, VA 23665. 
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1. INTRODUCTION 

The use of unstructured triangular meshes in two-dimensions and tetrahedral meshes in 
three dimensions is becoming more widely accepted as a viable means for computing compres- 
sible flows about complex geometries. Various on-going research efforts in this area can be 
found in the literature [1,2,3]. The main advantage of unstructured meshes is that they provide 
a flexible manner for discretizing topologically complex domains without recourse to large 
scale domain decomposition, which is necessary when structured meshes are employed [4]. 
They also provide a natural setting for the use of adaptive meshing. In fact, previous efforts of 
adaptive meshing on structured quadrilateral meshes have demonstrated the need for an 
unstructured-data approach [5]. 

On the other hand, unstructured mesh solvers have most often fallen short of their struc- 
tured mesh counteiparts in terms of efficiency and accuracy. The efficiency of these solvers is 
hindered by the use of indirect addressing required by random data-sets, and by generally less 
efficient solution algorithms. The problem due to indirect addressing can be somewhat allevi- 
ated by the use of hardware gather-scatter functions available on an increasing number of 
present-day supercomputers. Recent work by the author [6,7] has also resulted in an unstruc- 
tured multigrid algorithm, which has been shown to produce convergence rates comparable to 
those obtained with current structured mesh solvers. It would thus seem that unstructured mesh 
solvers are becoming increasingly competitive with structured mesh solvers. 

However, the accuracy of unstructured mesh solutions has seldom been able to match that 
delivered by structured mesh solvers. This issue has in fact been a subject of recent debate. 
In a recent report, Roe [8] argues that, in general, a piecewise linear discretization of the Euler 
equations on an arbitrary unstructured mesh results in locally first-order accuracy. Similarly, 
finite-element proofs of accuracy are most often limited to elliptic equations, and contain con- 
straints on the geometry of the elements. On the other hand, there is much practical evidence 
indicating that global second-order accuracy can be achieved on unstructured meshes for the 
Euler equations [9,10]. While much of the accuracy debate for the Euler equations centers 
around the discretization of the convective terms on arbitrary meshes, it is found in practice 
that the solution accuracy is most significantly affected by: 

a) The "quality" of the mesh (smoothness, number of connections at a point) 

b) The construction of the artificial dissipation terms 

c) The formulation of the inner and outer boundary conditions 


It is the intent of this work to demonstrate how careful attention to these areas can pro- 
duce solutions of equivalent accuracy with regards to current structured mesh solvers. An adap- 
tive meshing procedure is also described, which produces smoothly varying high resolution 
meshes. When combined with the unstructured multigrid approach, this method produces 
efficient and accurate solutions of flows around complex configurations. 


2. DISCRETIZATION OF THE GOVERNING EQUATIONS 

The variables to be determined are the pressure, density, Cartesian velocity components, 
total energy and total enthalpy denoted by p, p, u, v, E, and H, respectively. Since for a perfect 
gas we have 
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where y is the ratio of specific heats, we need only solve for the four variables p, p u, pv, and 
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P E. 

These values are determined by solving the Euler equations, which in integral form read: 

■Jjr , dxdy + j (fdy - gdx ) = 0 (2) 

where Q is a fixed area with boundary 3Q, x and y are Cartesian coordinates, and 
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The w variables are stored at the vertices of each triangle. The control volume for a vertex i is 
defined as the union of all triangles having a vertex at i, as shown in Figure 1. The boundary 
flux integral in equation (2) is approximated by first calculating the values of the fluxes f and g 
at the nodes on the outer boundary of this control volume. These can then be integrated around 
the control volume boundary by assuming that on each edge, the value of the flux can be taken 
as the average of the two values on either end of the edge. This corresponds to a finite- 
volume formulation using the trapezoidal rule of integration. It can also be shown to be 
equivalent to a piecewise linear Galerkin approximation, with a lumped mass matrix [1], 

Additional dissipative terms are required. These are formulated as a blend of an undi- 
vided Laplacian operator in the flow variables, and an undivided biharmonic operator. The 
biharmonic dissipation is formally third-order accurate on regular meshes and is applied 
throughout the flow-field. The Laplacian dissipation is formally first-order accurate on regular 
meshes and is applied only in the vicinity of a shock, where it is needed to prevent oscilla- 
tions. This adaptive control of the Laplacian dissipation is achieved by multiplying it by a 
Laplacian of the pressure, which is of order unity near a shock, and remains small elsewhere. 
This formulation of the dissipation represents an extension of the blended second and fourth 
differences originally proposed by Jameson, Schmidt, and Turkel [11]. While better shock 
resolution can be obtained with various upwind based methods [3,12], the present formulation 
is generally better suited for explicit schemes and less expensive to compute. This, in turn, 
suggests that increased accuracy may otherwise be obtained by adaptive refinement of the mesh 
in critical areas of the flow-field. 

3. INTEGRATION TO A STEADY-STATE 

Discretization of the Euler equations in space transforms the governing equations into a 
set of coupled ordinary-differential equations which must be integrated in time to obtain the 
steady-state solution. Thus, equation (2) becomes the set 

dW: 

Si-^- + [Q(w i )-D( Wi )] = 0, i= 1,2,3,... (4) 

where Si is the area of the control volume i, and is independent of time. The convective opera- 
tor Q(w) represents the discrete approximation to the flux integral in (2), and the dissipative 
operator D(w) represents the artificial dissipation terms. These equations are integrated in time 
using a fully explicit 5-stage hybrid time-stepping scheme, where the operator Q(w) is 
evaluated at each stage in the time step, and the operator D(w) is only evaluated in the first 
two stages, and then frozen at that value. Thus we advance in time as 

w (°) _ y/l 
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where w n and w n+1 are the values at the beginning and the end of the nth time step. The stan- 
dard values of the coefficients are 


a, = 1/4 a 2 = 1/6 a 3 = 3/8 a 4 = 1/2 a 5 = 1 

This scheme represents a particular case of a large class of hybrid time-stepping schemes, 
which has been specifically designed to produce strong damping characteristics of high- 
frequency error modes. This is a desirable property when a multigrid algorithm is to be 
employed. 

The unstructured multigrid algorithm, which has previously been described in detail 
[6,7,9], is used to accelerate convergence of the fine mesh solution to a steady state. This is 
achieved by repeatedly time-stepping on a sequence of coarse and fine unstructured meshes. 
The various meshes of the sequence are assumed to be independent of one another, and the 
pattern of transferring variables back and forth between the various meshes is determined in a 
preprocessing stage, where an efficient tree-search algorithm is employed. Convergence to a 
steady-state is also accelerated by using the maximum permissible time step at each point in 
the flow-field, as determined by local stability analysis, by the use of enthalpy damping [11], 
and implicit residual averaging [9]. 


4. FORMULATION OF THE DISSIPATION 

The artificial dissipation is constructed as a blend of undivided Laplacian and biharmonic 
operators. For subcritical flow, only biharmonic dissipation is required. This operator is con- 
structed by first forming an approximation to a Laplacian in the flow variables at each node 
according to : 

n n 

dWi = - Wi ] or dwi = Y, w k ~ nw { (6) 

*=1 Jt=i 


w represents the flow variables as given in equation (3), with the exception of the last term, 
where w is taken as the total enthalpy rather than the total energy, in order to preserve a con- 
stant total enthalpy in the steady-state solution, k represents the neighbors of node /, the total 
number of which is n, as shown in Figure 1. On a cartesian grid, this corresponds to the stan- 
dard five-point finite difference formula. The biharmonic dissipation is then formed by repeat- 
ing the above step: 


where A,- is given by 
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and e is an empirically determined coefficient. As in the structured mesh counterpart, this dissi- 
pation is independent of the time step. The two-dimensional Euler equations have four eigen- 
values, i. e. «, u, u+c, u-c, u being the local fluid velocity, and c the local speed of sound. The 
above construction of A t represents the maximum eigenvalue Iwl+c integrated around the boun- 
dary of the control volume for node i. This is consistent with the finite-volume formulation and 
also provides an isotropic value for scaling the dissipation with the maximum local eigenvalue 
about arbitrarily shaped control volumes. Furthermore, by averaging A t with A k , the dissipation 
is conservative in the w quantities. This can be seen by noting that for each contribution to the 
summation for node i in equations (6) and (7), an equivalent but opposite term will appear in a 
similar summation for node k. The Laplacian dissipation required for supercritical flows is con- 
structed in an analogous manner. 

5. BOUNDARY CONDITIONS 

The specific treatment of both convective terms and dissipative terms at both inner and 
outer boundaries was found to have a significant effect on the solution accuracy. For a vertex 
scheme, values are computed on the boundaries. The boundary control volumes resemble one 
half of a regular inner point control volume, as shown in Figure 2. 

For a vertex discretization of the Euler equations, wall boundary conditions may be 
applied in several manners. The weak form, consists of simply setting the normal components 
of the f and g fluxes in equations (3) to zero, except for the pressure terms in the momentum 
equations. In this formulation however, the velocities calculated on the boundary are not neces- 
sarily tangential to the wall. A stronger formulation consists of explicitly setting these veloci- 
ties to be tangential to the wall. However, for a discretized representation of the body, the 
tangent at a point is ill defined. In practice, a tangent can be defined as the line connecting the 
two neighboring boundary points. For curved boundaries on unequally spaced meshes, and for 
boundaries of varying curvature, this definition of the tangent will yield non-zero fluxes 
through the faces of the wall. Hence, in addition to setting the velocities tangential to the wall 
as defined above, the fluxes through the wall are also set to zero. In general, the resulting flux 
balance for the two momentum equations at the wall will not be satisfied. However, it is the 
combination of the flux balance and the velocity corrections at the boundary which sum to zero 
at convergence. In itself, this formulation of the inner boundary condition has only a moderate 
effect on the solution accuracy. However, the level of dissipation needed to guarantee conver- 
gence with this formulation is an order of magnitude smaller than with the weak formulation. 
It has been found that, as the level of dissipation is lowered, odd-even decoupling of the solu- 
tion does not occur. Rather, as the convergence degrades, the solution begins to exhibit low- 
frequency oscillations. Since high-frequency oscillations may be observed on triangular meshes 
when the weak formulation of the boundary condition is employed, and since they also occur 
when the strong formulation is employed on quadrilateral meshes, it appears that the combina- 
tion of the strong formulation and the increased connectivity of the triangular meshes inhibit 
odd-even decoupling of the solution. In fact, when decoupling is observed with the weak for- 
mulation of the boundary condition on triangular meshes, the wall velocity vectors are seen to 
oscillate about the tangent direction. Furthermore, whereas an odd-even oscillatory mode 
which exactly satisfies the discretized equations in the interior domain can be found for cell- 
centered and vertex discretizations on quadrilaterals, and for cell-centered discretizations on tri- 
angles [9], this is not possible for vertex discretizations on triangles. On a regular triangular 
mesh, as shown in Figure 3, only a three mode oscillation pattern can be found for a vertex 
scheme. In practice, the value of the e coefficient of equation (7) used in this work is one third 
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of the minimum value of the coefficient which would prevent high frequency oscillations when 
the weak formulation was employed. This represents a compromise between solution accuracy 
and speed of convergence. 

Treatment of the artificial dissipation terms at the inner boundaiy can also have a 
significant effect on the solution accuracy. Large amounts of numerical entropy may be gen- 
erated in regions of high gradients such as in the rapid expansion on the upper surface of an 
airfoil at high angles of attack. This false entropy contaminates the solution as it is convected 
downstream. In particular, when the upper and lower surface flows over an airfoil meet at the 
trailing edge, large oscillations may occur if one of the flows has an appreciable level of 
numerical entropy. To help minimize the spurious production of entropy by the boundary dissi- 
pation terms, it is convenient to divide the flow gradients into normal and tangential com- 
ponents at the boundary. For example, this is often achieved with body-fitted quadrilateral 
meshes where the mesh lines are perpendicular to the body. However, for a triangular mesh, it 
is impossible to have all mesh lines or edges hitting the boundary either tangentially or nor- 
mally. Thus assumptions concerning the local flow gradients at the boundary must be made. In 
the absence of body curvature, the boundary conditions for the Euler equations imply that the 
normal derivatives of all quantities at the wall vanish. Hence, it has been assumed that the nor- 
mal gradient is zero, when computing the artificial dissipation at the wall. Thus, for a mesh 
edge which intersects the wall at some prescribed angle, the difference is replaced by 
w\-w h where w, is the value at the wall, w k represents the value at the other end of the edge, 
and w' k is the value of the point on the wall obtained by projecting node k down onto the 
boundary, as shown in Figure 2. The specific value at this point is obtained by linear interpo- 
lation between the two points on the wall at either end of this interval. It can be shown that 
this construction is equivalent to that obtained by omitting the dissipation along all edges 
touching the boundary, and multiplying the tangential boundary dissipation values by a con- 
stant determined by the local grid topology. 

Global flow characteristics such as lift and drag coefficients for an airfoil can also be 
affected by the treatment of the outer boundary. Characteristic boundary conditions have been 
employed at the far-field boundary. The far-field boundary points of the unstructured meshes 
are not updated by the time-stepping procedure employed for interior points. Rather, the values 
at these points are updated by assuming the flow to be locally one-dimensional, and calculating 
the values at the next time step according to the ingoing and outgoing Rieman invariants for 
this simplified flow. For the outgoing Rieman invariant, values are extrapolated from neighbor- 
ing interior points. For the incoming Rieman invariant, flow values are taken as the superposi- 
tion of the freestream flow with the disturbance flow-field set up by a vortex centered at the 
quarter chord point of the airfoil, producing total circulation equal to that of the lifting airfoil 
[13]. Two other conditions are required for the full two-dimensional Euler flow. These are 
given by the values of tangential velocity and entropy, which are extrapolated from from inte- 
rior values for an outflow boundary, and set equal to to freestream values, augmented by the 
vortex flow-field, for an inflow boundary. Dissipation at the outer boundary is treated in the 
same manner as at the inner boundary, i.e. all normal gradients are assumed to vanish when 
computing the dissipation terms. This is equivalent to the treatment advocated for cell-centered 
quadrilateral Euler solvers. 

6. RESULTS ON SIMPLE MESHES 

The first case considered is a NACA0012 airfoil with a freestream Mach number of 0.8, 
and an incidence of 1.25°. The mesh employed for this case consists of a regular triangular grid 
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with 4224 nodes, as shown in Figure 4. The positioning of the nodes is derived from a 128 x 
32 quadrilateral O-mesh. This mesh represents the finest of a sequence of four meshes used in 
the multigrid algorithm to accelerate convergence. The computed surface pressure distribution 
is shown in Figure 5. The upper and lower shocks arc both well defined, and no oscillations 
are visible. The lift and drag coefficients compare well with other inviscid solutions on 
equivalent and finer quadrilateral meshes [14]. The surface entropy, which is shown in Figure 
6, is calculated as 



( 9 ) 


The entropy values are nearly zero everywhere ahead of the strong shock, with a maximum 
value near the leading edge of 0.3%. Behind the upper surface shock, the entropy assumes a 
constant value as it is convected downstream. The multigrid convergence history given in Fig- 
ure 7, indicates a reduction of the RMS average of the density residuals of nearly 6 orders of 
magnitude over 200 cycles. 

The second case is that of a Korn airfoil at design conditions. This airfoil is designed to 
produce shock-free supercritical flow at these conditions. Since such flows are known to be 
isolated points (i.e. any change in the geometry or freestream conditions produces a shock), the 
ability to capture such flows is a good indication of the accuracy of the scheme. The mesh 
employed corresponds to a 128 x 32 O-mesh and is depicted in Figure 8. From the computed 
pressure distribution in Figure 9, the compression of the flow through the transonic regime can 
be seen to be gradual. The low drag value (Cd = 0.0002) is an indication of the isentropic 
nature of the flow, as well as the accuracy of the computed solution. 

The third case involves flow about a Karman-Trefftz airfoil and flap. An exact 
incompressible potential solution exists for this case which may be used for comparison pur- 
poses. The unstructured mesh about this configuration is shown in Figure 10. In Figure 11, the 
calculated surface pressure distribution is compared with the exact solution. This case was run 
for a freestream Mach number of 0.125. Hence, effects due to compressibility are small (<1% 
according to the Prandtl Glauert rule). Good agreement between the two solutions is observed 
with respect to the surface pressure plots, and the calculated force coefficients, which are given 
in Table 1. 

7. ADAPTIVE MESHING STRATEGY 


For flows about more complicated geometries, the location of regions where rapid varia- 
tions in the flow occur are seldom known at the outset. Thus, an adaptive meshing strategy is 
necessary to ensure accurate solutions of such flows, as well as to maintain the efficiency of 
the solver. 

The use of adaptive meshing, in conjunction with the unstructured multigrid algorithm 
has previously been described in [7]. The method first computes the flow on an existing coarse 
mesh. A new finer mesh is then constructed by generating a new distribution of grid points as 
determined by the computed flow-field on the coarse grid. These points are then joined 
together to form an unstructured triangular mesh, using the Delaunay triangulation procedure. 
Delaunay triangulation is a geometrical construct which represents a unique triangulation of a 
given set of points in a plane. The resulting mesh is then smoothed by repositioning the mesh 
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points as detennined by a Laplacian filtering operation. The patterns for transferring variables 
back and forth between this new mesh and the previous coarser mesh are determined prior to 
the solution of the flow-field using a tree-search algorithm. This approach does not assume 
any relation between the two meshes, thus offering great flexibility in configuring the finer 
mesh for optimum accuracy. The flow variables are then transferred to the new mesh and 
time-stepping proceeds, employing the previous mesh as a coarse mesh in the multigrid 
sequence. This procedure may be repeated, each time generating a new finer mesh, until the 
desired level of accuracy is achieved. At present, new mesh point distributions are generated 
by combining all previous mesh points with a set of newly generated points. The distribution 
of these new mesh points is determined according to a criterion developed by Dannenhoffer 
and Baron [15]. The first undivided difference of the density along a mesh edge is examined. 
When this value is larger than some fraction of the RMS average of all density differences 
across the mesh, a new point is added midway along the edge. A second pass is then per- 
formed which splits the remaining edges of each triangular element bordering on a previously 
split edge, thus ensuring an isotropic refinement in these regions. A third pass of this type is 
also performed on the boundary elements. The mesh smoothing operation described above dis- 
places the mesh points and thus, in general, the coarse grid points will not coincide with any of 
the fine grid points. Previous efforts at directional refinement by refining only along edges 
have met with limited success, due to the nature of Delaunay triangulation, which tends to 
minimize cell aspect ratios. Other efforts in the literature at directional refinement using alter- 
nate triangulation procedures have been more successful [3]. 

8. RESULTS EMPLOYING ADAPTIVE MESHES 

A subcritical and a supercritical case about a high-lift three-element airfoil configuration 
were computed using the adaptive meshing strategy. These solutions were compared with a full 
potential flow solver made available by The Boeing Company [16]. This method is a finite- 
element non-linear potential flow solver which is formally first-order accurate in supersonic 
zones, and second-order accurate in subsonic zones. The method operates on a structured qua- 
drilateral mesh, which is generated from a streamline pattern obtained with the aid of a prelim- 
inary panel method solution. 

In the subcritical case, the freestream Mach number is 0.2, and the incidence is 8°. The 
final mesh employed in this calculation is shown in Figure 12. It contains a total of 10745 
points, of which 432 are on the airfoil surfaces. This represents the sixth mesh level of the 
multigrid sequence in which the first two levels are generated globally as described in [6,9], 
and the four finest levels are generated adaptively. A globally generated mesh of this resolution 
would have required 10 to 20 times more points, and thus would be prohibitively expensive. 
The computed surface pressure distributions are compared with those obtained by the full 
potential solution in Figure 13. This case exhibits an extreme double suction peak near the 
leading edge of the main airfoil. Good agreement is observed between the Euler and potential 
solutions, including the predicted magnitude of the suction peaks. For subcritical cases, the 
flow should be isentropic, and thus, the plot of the surface entropy in Figure 14, can be inter- 
preted as a measure of the error in the solution. This entropy is nearly zero everywhere, 
except for a region near the leading edge of the main airfoil, where it reaches a peak of 2 %. 
However, throughout most of this region, the entropy is less than 1%. The maximum entropy 
level could be reduced to less than 1% by employing an additional mesh level. However, the 
effect on the surface pressure distribution was found to be minimal. The convergence rate for 
this case is depicted in Figure 15, where a residual reduction of seven orders of magnitude on 
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the finest mesh was achieved over 300 multigrid cycles, employing all six meshes in the mul- 
tigrid sequence. A plot of the Mach contours in the flow-field is also given in Figure 16. 

For the supercritical case, the Mach number is increased to 0.25, while the incidence 
remains 8°. At these conditions, the flow becomes supersonic as it expands around the leading 
edge of the main airfoil, and a shock is formed. Seven mesh levels were used to compute this 
solution. The finest mesh is illustrated in Figure 17, where the highest level of refinement is 
seen to occur mainly near the leading edge of the main airfoil. The total number of grid points 
is 12830, of which 539 are on the airfoil surfaces. A convergence of five orders of magnitude 
over 300 multigrid cycles was achieved for this case, employing seven meshes in the multigrid 
sequence. The computed surface pressure distribution is shown in Figure 18, where it is com- 
pared to the potential flow solution. Figures 20 and 21 also depict the Mach number contours 
in the flow-field. Although the extent of the shock is small, the Mach number ahead of it is 
1.53 in the Euler solution, or 1.60 in the potential flow solution. In either case, this represents 
a strong shock, which cannot be accurately modeled by the potential flow assumption. 
Although there are some discrepancies between the two solutions on the main airfoil, the larg- 
est differences occur on the upper surfaces of the vane and flap. This is surprising, since the 
entirety of the flow which impinges upon the vane and the flap originates from the lower sur- 
face of the main airfoil, as it is channeled through the gaps. It appears that the shear layer pro- 
duced by the shock on the upper surface of the main airfoil (see Figure 20), which is not 
present in the potential solution, has a global effect on the entire flow-field. On the other hand, 
the numerical diffusion of this shear layer, which enables it to impinge slightly upon the vane 
and flap, may contribute to these discrepancies. However, from Figure 19, the surface entropy 
is seen to be much less than 1% everywhere except in the shear layer behind the shock, on the 
main airfoil. - Furthermore, the same characteristic discrepancies between the two solutions 
were noted when alternate meshes were employed for the Euler solution. Increased mesh reso- 
lution in the region of the shear layer would be desirable to further limit the diffusion of 
entropy. 

9. CONCLUSION 

An accurate and efficient method for solving the Euler equations about complex 
geometries has been demonstrated. The method requires approximately 0.12xl0 -3 
secs/node/mg-cycle on the CRAY-2 supercomputer. Thus, the results about simple geometries 
presented in this work required 1 to 2 minutes of CPU time, whereas the three-element airfoil 
calculations were obtained in about 6 to 8 minutes. The present code can be expected to run 
substantially faster on a CRAY-XMP class machine, due to the more efficient gather-scatter 
routines available. Good agreement with other available inviscid solutions has been demon- 
strated. For supercritical flow over the three-element high-lift configuration, discrepancies 
between the present Euler solution and the full potential solution of [16] exist, the origin of 
which appear to be related to the presence of a shear layer behind the shock, which cannot be 
captured by a full potential method. Future work is aimed at further optimizing the efficiency 
of the method by implementing a directional refinement strategy. The inclusion of viscous 
effects by solving the Reynolds-averaged Navier-Stokes equations is also planned. 
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Figure 1 

Control volume for vertex triangle scheme 




Figure 3 

Three-Mode Decoupling Possible on a Regular Triangular Mesh 
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Figure 4 

Partial View of Mesh about NACA0012 airfoil 
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Efemt#: 1 Nnode: 4460 Ncyc: 200 
cl : 0.3543 cd : 0.0224 cm :-0.0373 


Figure 5 

Computed Surface Pressure Distribution for NACA0012 Airfoil 
Mach = 0.8, Incidence = 1.25* 
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Elemt#: 1 Nnode: 4460 Ncyc: 200 
cl : 0.3543 cd : 0.0224 cm :-0.0373 


Figure 6 

Computed Surface Entropy for NACA0012 Airfoil 
Mach = 0.8, Incidence = 1.25° 
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Residl : 0.12e+01 Resid2 : 0.37e-05 Rate : 0.9381 


Figure 7 

Convergence History for the supercritical NACA0012 case 
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Figure 8 

Partial View of Mesh about the Korn Airfoil 
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Elemt#: 1 Nnode: 4480 Ncyc: 400 
cl : 0.6302 cd : 0.0002 cm :-0.1463 


Figure 9 

Computed Surface Pressure Distribution for the Korn Airfoil at Design Conditions 
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Main Airfoil 

Hap 

Total 


Cl 

Cd 

Cl 

Cd 

Cl 

Cd 

Nodal 

Scheme 

1.7005 

-0.0906 

0.3357 

0.0890 

2.0362 

-0.0016 

Analytical 

Incompressible 

Solution 

1.6915 

-0.0898 

0.3366 

0.0897 

2.0281 

-0.0001 


Table 1 

Comparison of Computed Force Coefficients for the Karman-Trefftz Configuration 





















Figure 12 

Illustration of the Adaptively Refined Mesh for the Subcritical Case of the 3-Element Airfoil Configuration 

10745 nodes, 20987 triangles 
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Figure 14 

Computed Surface Entropy for the Subcritical 3 -Element Airfoil Case 
Mach = 0.2, Incidence = 8° 
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Residl : 0.35e-h01 Resid2 : 0.42e— 06 Rate : 0,9481 


Figure 15 

Convergence Rate on Finest Mesh for the Subcritical 3-Element Airfoil Case 
as Measured by the RMS Average of the Density Residual throughout the Flowfield 
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Figure 16 

Flowfield Mach Contours for the Subcritical 3-Element Airfoil Case 
Mach = 0.2, Incidence = 8° 
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Figure 17 

Illustration of the Adaptively Refined Mesh for the Supercritical Case of the 3 -Element Airfoil Configuration 

12830 nodes, 25125 triangles 
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Figure 18 

Comparison of Computed Surface Pressure Distribution for the Euler Solution 
with that of the Full Potential Solution; Mach = 0.25, Incidence = 8° 
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Figure 20 

Flowfield Mach Contours for the Supercritical 3-Element Airfoil Case 
Mach = 0.25, Incidence = 8° 




Figure 21 

Mach Contours near the Leading Edge of the Main Airfoil of the 3-Element Configuration 
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